#ifndef __ARRAY1D_H__
  #define __ARRAY1D_H__

#include <string>
#include <sstream>
#include <climits>

#include "exceptions.h"
#include "macros.h"
#include "global_variables.h"
#include "useful_functions.h"

using namespace std;

//! One dimensional array

/*!
A simple array implementation.
For arrays larger than 1GB some changes may be necessary.

Implementation date: 2008

Usage:

\code
array1d<int> a(3); // with 3 elements

a(0) = 17;

a(1) = 21;

a(2) = -4;
\endcode
*/

template <class T>
class array1d {
private:
  T *content;
  unsigned int stream_position;
  unsigned int rows;

  void allocate();
  void release();
  void hard_copy(const array1d &a);

public:
  //! resets the array to an empty array and releases allocated memory
  void clear();

  //! initializes all elements of the array to init_value
  void initialize(const T &init_value);

  //! returns the position in which the next element will be inserted when operator<< is called
  unsigned int get_stream_position() const;

  //! sets the position in which the next element will be inserted when operator<< is called
  void set_stream_position(unsigned int position);

  //! returns the string representation of the object
  string to_string() const;

  //! returns the smallest element
  const T &get_minimum() const;

  //! returns the largest element
  const T &get_maximum() const;

  //! returns the sum of all elements
  T get_sum() const;

  //! returns the number of rows
  unsigned int get_row_count() const;

  //! reconstructs the array and deletes its content
  void reconstruct(unsigned int rows);

  //! resizes the array and keeps its content intact
  void resize(unsigned int rows, const T &init_value);

  //! sorts the array using quick sort
  void sort();

  //! searches for a value inside the array and returns its index or UINT_MAX if not found
  unsigned int search(const T &value) const;

  //! searches for a value inside the array using binary-search (the array must be already sorted) and returns its index or UINT_MAX if not found
  unsigned int binary_search(const T &value) const;

  //! negates and returns a copy of the array
  array1d negate() const;

  //! mathematical addition (may throw an exception)
  array1d operator+(const array1d &a) const;

  //! mathematical subtraction (may throw an exception)
  array1d operator-(const array1d &a) const;

  //! mathematical multiplication
  array1d operator*(const T &value) const;

  //! index operator (may throw an exception)
  T &operator()(unsigned int row) const;

  //! assignment operator
  array1d &operator=(const array1d &a);

  //! shift operator: inserts an element into the array at the stream_position (may throw an exception)
  array1d<T> &operator<<(const T &value);

  //! copy constructor
  array1d(const array1d &a);

  //! default constructor: an empty array, reconstruct or resize must be called to enlarge the array before it can be used
  array1d();

  //! constructor: does not initialize the array content
  explicit array1d(unsigned int rows);

  //! constructor: initializes the array content
  explicit array1d(unsigned int rows, const T &init_value);

  //! destructor (non-virtual: please do not inherit from this class)
  ~array1d();
}; // class array1d

//////////////////////////////////////////////////////////////////////////////

template <class T>
void array1d<T>::allocate() {
  release();
  content = new T[rows];
}

template <class T>
void array1d<T>::release() {
  delete[] content;
  content = 0;
}

template <class T>
void array1d<T>::hard_copy(const array1d &a) {
  this->stream_position = 0;
  this->rows = a.rows;

  allocate();

  #pragma omp parallel for if (rows > MINIMUM_SIZE_LIMIT_OMP)
  for (unsigned int i = 0; i < rows; ++i) // do not use memcpy (it may cause problems with non-standard datatypes)
    this->content[i] = a.content[i];
}

template <class T>
void array1d<T>::clear() {
  release();
  stream_position = 0;
  rows = 0;
}

template <class T>
void array1d<T>::initialize(const T &init_value) {
  #pragma omp parallel for if (rows > MINIMUM_SIZE_LIMIT_OMP)
  for (unsigned int i = 0 ; i < rows; ++i)
    content[i] = init_value;
}

template <class T>
unsigned int array1d<T>::get_stream_position() const {
  return this->stream_position;
}

template <class T>
void array1d<T>::set_stream_position(unsigned int position) {
  this->stream_position = position;
}

template <class T>
string array1d<T>::to_string() const {
  ostringstream result;

  for (unsigned int i = 0 ; i < rows; ++i) {
    result << '[';
    result << i;
    result << "] = ";
    result << content[i];
    result << '\n';
  } // for i

  return result.str();
}

template <class T>
const T &array1d<T>::get_minimum() const {
  if (rows == 0)
    throw operation_not_allowed_exception("Array is empty!");

  return *min_element(content, content + rows);
}

template <class T>
const T &array1d<T>::get_maximum() const {
  if (rows == 0)
    throw operation_not_allowed_exception("Array is empty!");

  return *max_element(content, content + rows);
}

template <class T>
T array1d<T>::get_sum() const {
  if (rows == 0)
    return 0;

  T total = content[0];

  #pragma omp parallel for if (rows > MINIMUM_SIZE_LIMIT_OMP) reduction(+ : total)
  for (unsigned int i = 1; i < rows; ++i)
    total += content[i];

  return total;
}

template <class T>
unsigned int array1d<T>::get_row_count() const {
  return rows;
}

template <class T>
void array1d<T>::reconstruct(unsigned int rows) {
  this->rows = rows;

  allocate();

  set_stream_position(0); // makes sure that operator<< will begin from the index 0
}

template <class T>
void array1d<T>::resize(unsigned int rows, const T &init_value) {
  if (this->rows >= rows)
    return; // do not shrink the array now!

  T *tmp = new T[rows];

  #pragma omp parallel for if (rows > MINIMUM_SIZE_LIMIT_OMP)
  for (unsigned int i = 0 ; i < rows; ++i) {
    if (i < this->rows)
      tmp[i] = content[i];
    else
      tmp[i] = init_value;
  } // for i

  this->rows = rows;

  release();
  content = tmp;
}

template <class T>
void array1d<T>::sort() {
  qsort(content, rows, sizeof(T), comparison<T>);
}

template <class T>
unsigned int array1d<T>::search(const T &value) const {
  for (unsigned int i = 0 ; i < rows; ++i) {
    if (content[i] == value)
      return i;
  } // for i

  return UINT_MAX; // not found
}

template <class T>
unsigned int array1d<T>::binary_search(const T &value) const {
  const T *p = static_cast<const T*>(bsearch(&value, content, rows, sizeof(T), comparison<T>));
  if (p == 0)
    return UINT_MAX; // not found
  else
    return static_cast<unsigned int>(p - content);
}

template <class T>
array1d<T> array1d<T>::negate() const {
  return (*this) * -1;
}

template <class T>
array1d<T> array1d<T>::operator+(const array1d &a) const {
  if (a.rows != this->rows)
    throw operation_not_allowed_exception(global_empty_string);

  array1d result(*this);

  #pragma omp parallel for if (rows > MINIMUM_SIZE_LIMIT_OMP)
  for (unsigned int i = 0 ; i < rows; ++i)
    result.content[i] += a.content[i];

  return result;
}

template <class T>
array1d<T> array1d<T>::operator-(const array1d &a) const {
  return (*this + a.negate());
}

template <class T>
array1d<T> array1d<T>::operator*(const T &value) const {
  array1d result(*this);

  #pragma omp parallel for if (rows > MINIMUM_SIZE_LIMIT_OMP)
  for (unsigned int i = 0 ; i < rows; ++i)
    result.content[i] *= value;

  return result;
}

template <class T>
T &array1d<T>::operator()(unsigned int row) const {
  INDEX_CHECK(row, this->rows);

  return content[row];
}

template <class T>
array1d<T> &array1d<T>::operator=(const array1d &a) {
  if (&a == this)
    return *this;

  hard_copy(a); // content will be released first
  return *this;
}

template <class T>
array1d<T> &array1d<T>::operator<<(const T &value) {
  INDEX_CHECK_MSG(stream_position, rows, "Out of range! enlarge the array first.");

  content[stream_position] = value;
  ++stream_position;

  return *this;
}

template <class T>
array1d<T>::array1d(const array1d &a) : content(0), stream_position(0), rows(a.rows) {
  hard_copy(a);
}

template <class T>
array1d<T>::array1d() : content(0), stream_position(0), rows(0) {
  // nothing to do
}

template <class T>
array1d<T>::array1d(unsigned int rows) : content(0), stream_position(0), rows(rows) {
  allocate();
}

template <class T>
array1d<T>::array1d(unsigned int rows, const T &init_value) : content(0), stream_position(0), rows(rows) {
  allocate();
  initialize(init_value);
}

template <class T>
array1d<T>::~array1d() {
  release();
}

#endif
